A novel FLEX supplemented QMC approach to the Hubbard model 
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This paper introduces a novel ansatz-based technique for solution of the Hubbard model over two 
length scales. Short range correlations are treated exactly using a dynamical cluster approximation 
QMC simulation, while longer-length-scale physics requiring larger cluster sizes is incorporated 
through the introduction of the fluctuation exchange (FLEX) approximation. The properties of the 
resulting hybrid scheme are examined, and the description of local moment formation is compared 
to exact results in ID. The effects of electron-electron coupling and electron doping on the shape of 
the Fermi-surface are demonstrated in 2D. Causality is examined in both ID and 2D. We find that 
the scheme is successful if QMC clusters of Nc > 4 are used (with sufficiently high temperatures in 
ID), however very small QMC clusters of Nc = 1 lead to acausal results. 

PACS numbers: 71.10.Hf 



I. INTRODUCTION 

Despite its transparent meaning and simple form, the 
Hubbard model is thought to contain much of the funda- 
mental physics of correlated electron systems 0. Since 
its inception, the Hubbard model has been successfully 
applied to the description of the Mott-Hubbard metal- 
insulator transition, and is thought to contain a minimal 
description of certain transition metal oxides, such as the 
high-Tc cuprates and the low-Tc ruthcnatcs. While lim- 
iting aspects of the model are easily understood, the ex- 
act ground state of the one-band Hubbard model remains 
unknown for all cases other than the Hubbard chain 0|. 

Examination of the infinite-dimensional limit of the 
Hubbard model using dynamical mean-field theory 
(DMFT) has lead to a greater understanding of inter- 
mediate coupling phenomena such as the Mott-Hubbard 
metal-insulator transition DMFT maps the full lat- 
tice problem onto a single impurity Anderson model, 
which may be solved using various numerically exact 
techniques such as quantum Monte-Carlo (QMC) [4J and 
numerical renormalization group (NRG) 5]. In recent 
work, the dynamical mean-field theory was extended to 
study correlated electron systems in 2D, resulting in an 
approach known as the dynamical cluster approximation 
(DCA) 0. DCA is a fully causal method, which system- 
atically restores non-local corrections to the dynamical 
mean-field theory. In previous work, Jarrell et al. com- 
puted the phase diagram of the 2D Hubbard model for 
small clusters of Nc = 4 using QMC 0. Unfortunately, 
the QMC solution of large clusters is prohibitively expen- 
sive in terms of super-computer time, while alternative 
approaches, such as the perturbation-theory-based FLEX 
approximation can be used to solve large DCA clusters, 
but only in certain limits. FLEX is physically intu- 
itive, concentrating on scattering from important mecha- 
nisms (spin fluctuation, density fluctuation and particle- 
particle pairs), and in the past, FLEX has been shown 
to reproduce long-length-scale physics of the Hubbard 
model such as the Mermin- Wagner theorem 



In this paper we present an ansatz for combining 
long length scale information from the FLEX approxi- 
mation with the QMC solution of the Hubbard model. 
The two techniques are complimentary, since QMC pre- 
dicts the correct short length scale physics, while FLEX 
shows appropriate long length-scale behavior. In section 
UTI we introduce the DCA and the ansatz for the self- 
energy, and provide technical details of the modified self- 
consistent procedure for the incorporation of two-length- 
scale physics. We examine issues of causality and appli- 
cability in section II I II In section IIVI we examine local 
moment formation and the momentum dependent occu- 
pation nk- Comparisons are made between the hybrid 
technique introduced here, exact results for the ID model 
and the more traditional FLEX technique. We calculate 
the Fermi surface of the 2D model in section Finally, 
we discuss the applicability and outlook for the new tech- 
nique in section lv"Tl 



II. FORMALISM 

The DCA extends the DMFT, which assumes that the 
self-energy is constant across the Brillouin Zone (BZ), 
by introducing a limited momentum dependence corre- 
sponding to short range spatial fluctuations. This is 
achieved by breaking up the Brillouin zone into a series 
of subzones, within which it is assumed that the self- 
energy is constant. This allows a momentum integrated 
(or coarse grained) Green function to be defined in an 
analogous manner to the definition of DMFT: 
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The Ki are defined as in reference |6| as the mean mo- 
mentum for a coarse grained cell. Via an inverse Dyson 
equation, this leads to a "bare" Green function corre- 
sponding to the host of a cluster impurity problem, 



E(K, nu n ) = ^(K, ito n ) - G- X (K, iu n ) . 



(2) 
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which may be solved with a variety of methods. The 
approximation has two well defined limiting cases. For 
cluster sizes of 1, the DC A maps onto the DMFT. For 
infinite subzones, the formalism is exact. For a finite 
number of subzones, results are expected to be closer to 
the behavior of the exact infinite lattice than conven- 
tional finite size techniques, since an infinite number of 
k-states is considered. The coarse graining over k-states 
maps the lattice model to be studied (in this case the 
Hubbard model) onto a periodic cluster in real space. 

In order to solve the resulting problem, approximations 
are employed. Two common cluster solvers are the FLEX 
approximation |9( , and the QMC method of Hirsch and 

I 



E^)(K,«j n ) = Eg^(K,ta; B ) 



The superscripts indicate the size of the cluster. Bars 
over the self energy indicate that a secondary coarse 
graining of the lattice self energy (linear interpolated self 
energy) from the other cluster size has been performed, 
i.e. 

E(K',iu, n ) = J2 S (M"V) ( 5 ) 

where E(k, iu) n ) is the linear interpolation (lattice self- 
energy) of the small cluster (an equivalent expression ex- 
ists to coarse grain from the large to the small cluster). 

The physical content of the ansatz can be clarified by 
examining the diagrammatic expansion in figure^ Since 
the QMC self energy is non-perturbative, it contains com- 
plete contributions from all orders in U for a cluster of 
size Nc (1st term). The second term removes all FLEX 
contributions for cluster size Nc- These are then replaced 
with the FLEX contributions for the larger cluster size 
N' c . This approach is reasonable if the set of FLEX dia- 
grams have a stronger momentum dependence than the 
remaining diagrams. The complete FLEX with particle- 
particle scattering obeys the Mermin-Wagner theorem 
in 2D, indicating a significant momentum dependence 
which partially justifies this assumption. 

As there are two length scales represented in the 
ansatz, two coarse grainings and cluster exclusions are 
performed consecutively (one for each length scale). The 
cluster excluded Green function for the small cluster is 
then used as input to the QMC, and the coarse grained 
Green functions are used as input to the FLEX. This 
leads to an iterative procedure which is demonstrated in 
figure |21 There are some additional methods that can be 
employed to aid convergence. The FLEX self-energy is 
strongly damped between iterations to avoid the FLEX 



Fye [H|T3- In principle, the QMC method gives the 
exact solution of the cluster, but the numerical cost is 
high, and only small clusters can be investigated. De- 
tails of the QMC and FLEX methods may be found in 
references and H respectively. 

In the hybrid technique proposed here, we define clus- 
ters of two sizes, one large of size N' c solved with the 
FLEX, and one small of size Nc solved with QMC. For 
the small cluster, coarse grained points are denoted K, 
and for the large cluster, they are denoted as K'. In or- 
der to couple the self energies from these methods, we 
define an ansatz 1131: 



4S(K, + fi(K, (3) 

E^ x (K',io; n ) + E^(K',^ n ) (4) 
I 

instability, and many FLEX iterations ~ 100 are carried 
out for each QMC step (the ansatz is recomputed on 
each sub-iteration of the FLEX). The QMC step is also 
damped to a lesser degree. 



III. RANGE OF APPLICABILITY 

Causality is essential for any predictive theory. Causal- 
ity violations in single-particle quantities can lead to neg- 
ative parts of the single-particle spectrum and density of 
states and violations of sum rules. Causality violations in 
two-particle quantities can lead to erroneous predictions 
for phase transitions. This section examines how causal- 
ity problems can arise, the expected regions of applica- 
bility for the new approximation, and tactics that may 
be employed to avoid numerical instabilities. Causality 
is reflected in the analytic properties of the self energy 
and Green functions. For example, the retarded self en- 
ergy E(k, ui) is analytic in the upper half of the complex 
frequency plane and — i±mE(k, ui) > 0. Since there is 
a subtraction in the ansatz for the self-energy, causal- 
ity problems might be anticipated, since the imaginary 
part of the real-frequency self-energy can become posi- 
tive. Here, we search for causality violations in the Mat- 
subara self energy. It is related to E(k, ui) through 

f -ilmE(k,w) 

E(k,zw„)= / dui^- - (6) 

J iui n - u) 

Consequently, — w n ImE(k, iui n ) > as a consequence of 
causality. We employ violations of this inequality as a 
sufficient but not necessary indication of the causality 
violation of our formalism. 
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FIG. 1: Feynman diagrams for the Hubbard model to 4th or- 
der. Shaded diagrams are computed for a large cluster, while 
the remaining diagrams are computed for the smaller clus- 
ter. In our scheme, an infinite number of FLEX diagrams is 
considered, and the remaining diagrams are calculated using 
QMC, which is a non-perturbative method. 
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FIG. 2: Flow chart showing the self-consistent procedure for 
the ansatz. Iteration is initialized using second order per- 
turbation theory (SOPT). The flow then continues with the 
calculations for both cluster sizes, which are carried out con- 
secutively. Once convergence is reached, analysis is carried 
out on the results from the large cluster to compute quantities 
of interest, such as the local moment and the Fermi-surface. 
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FIG. 3: Self-energy ansatz by contribution, including the 
QMC part for the small cluster (filled circles with solid lines), 
the FLEX part for the small cluster (filled squares with solid 
lines), the FLEX part for the large cluster (open squares with 
dashed lines) ,and the combination of these into the ansatz for 
the large cluster (open circles with dashed lines). U = O.SW, 
T = 0.08, N c = 4, N'o = 32, D = 1, n = 0.5. For each of the 
self-energy contributions, there are a number of curves repre- 
senting the different contributions at each momentum point 
as the Brillouin zone is crossed. Causality is preserved for 
this value of U. FLEX overestimates the self-energy, while 
QMC has the correct order of magnitude, but has insufficient 
details of the momentum dependence. The ansatz corrects 
the momentum dependence, while keeping the value of the 
self-energy at the appropriate order of magnitude. 



In figure |21 we examine the self-energy in ID for 
N c = 4, N' c = 32, U = 0.8W and T = 0.08 at half-filling. 
All components of the ansatz are shown: The QMC self- 
energy for the smaller cluster (filled circles), the FLEX 
self-energy for the smaller cluster (filled squares) and the 
FLEX self-energy for the larger cluster (open squares). 
The result of the ansatz is also shown for the larger clus- 
ter only (open circles). The combination of ID lattice 
and intermediate coupling results in an extreme test case 
for the ansatz, and causality can be seen to hold for both 
large and small clusters. 

There are two ways in which causality may be violated. 
The first relates to the nature of the ansatz. In order to 
avoid overcounting of diagrams, FLEX diagrams on the 
small length scale are subtracted before new diagrams 
on a larger length scale are reinserted. Therefore, even 
when the FLEX method gives causal results, the ansatz 
may return a non-causal self energy. The FLEX parts 
of the ansatz are designed to return the correct momen- 
tum dependence to the self energy, and may be seen as a 
perturbation to the QMC segment of the self energy, 



S(k) = E. 



QMC 



AE(k) 



(7) 



In order to avoid causality problems resulting from the 
subtraction of FLEX diagrams, small-length-scale clus- 
ters should be as large as possible, with a minimum re- 
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FIG. 4: Matsubara axis self energy calculated using FLEX. 
The top panel shows the self-energy for U = W, D = 1 and 
the bottom panel for U = 0.51^, D — 2. In both cases, the 
result from large and small clusters are sufficiently close to 
one another to avoid causality problems in the ansatz. 

quirement that Nc > 1. The reason for this is clear 
for very large QMC clusters, since in the limit that 
Nc — > N c , both FLEX solvers return the same results, 
and A£(k) — > i.e. causality is guaranteed. For a good 
cancellation of the FLEX terms, the small and large clus- 
ter FLEX self-energies should both be on the large N 
scaling curve. We demonstrate the scaling behavior of 
the FLEX approximation in figure 0] (a pure FLEX cal- 
culation is used). The imaginary part of the self-energy 
from the FLEX approximation is shown. Both ID and 
2D calculations are carried out, and Nc chosen to be 
representative of cluster sizes used in this paper. For the 
2D case, the Nc = 4 self-energy doesn't sit directly on 
the scaling curve. However, it has the correct form to 
promote causality since the extremal behavior is similar 
to that of the large cluster. In general, for smaller U, 
less cluster points are needed to fall on the scaling curve 
(with a minimum of Nc = 4 for both ID and 2D sys- 
tems). For larger U > W, the minimum Nc required to 



see correct scaling behavior in the FLEX approximation 
is expected to increase. We suggest that the use of the 
ansatz should be limited to couplings that are no bigger 
than the bandwidth. 

If the self-energy contributions from the FLEX approx- 
imation lie far away from the scaling curve, then the 
magnitude of A£(k) can be as large as the QMC part. 
FLEX inherently overestimates the magnitude of the self- 
energy, so errors are expected to be amplified at strong 
coupling. For this reason, use of the DMFT (Nc = 1) as 
a small cluster solver is not advised unless the tempera- 
ture is high, or the coupling is weak. 

Causality problems may also originate from the FLEX 
instability. FLEX is constructed from a geometric se- 
ries, which imposes a condition on the susceptibilities: 
Xpp(0, 0), Xph(Q, 0) < 1. At very high temperatures, 
both susceptibilities are small, and the geometric con- 
dition is met. As temperature decreases, the suscepti- 
bilities grow. In ID and 2D, where the Mcrmin- Wagner 
theorem holds, the FLEX instability should only cause 
minor issues of numerical instability, since the suscepti- 
bilities are never expected to diverge (N.B. Some numer- 
ical effort is still required to avoid the instability). In 
3D, however, FLEX predicts a phase transition at mod- 
erate temperatures. Below that transition, the unbro- 
ken symmetry approximation is neither causal nor valid. 
Lower temperatures could only be accessed by extend- 
ing the approximation to include the anomalous Green 
functions associated with broken symmetry states. It is 
unlikely that the Stoner criterion Xph = 1 is meaningful 
in the ansatz scheme, so there will be regions of the pa- 
rameter space that cannot be reached in 3D. In order to 
investigate the 3D phase diagram with this scheme, it is 
essential that the FLEX instability occurs at lower tem- 
peratures than the true phase transition. Alternatively a 
large cluster solver such as 2nd order perturbation theory 
could to be applied. The phase transition would be mea- 
sured through the divergence of the relevant 2-particle 
susceptibility. To calculate this quantity, it would be 
necessary to introduce an additional ansatz for the irre- 
ducible vertex function. This would be quite involved, 
and is left for a future publication. 

Finally, we note that FLEX (as an extension of the 
T-matrix approximation) is exact for dilute systems, and 
the ansatz is therefore expected to be exact in the region 
of low doping. So long as reasonable temperatures are 
maintained, results should be applicable to a wide range 
of couplings and dopings. For dopings closer to half- 
filling, the results of testing show the approximation to 
be valid for couplings of up to the band width. 



IV. ID MODEL 

In this section, the applicability of the hybrid 
FLEX/QMC approach to the ID Hubbard model is dis- 
cussed. Since the exact ground state of the ID case is 
known from the Bethe ansatz solution Q , a quantitative 
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FIG. 5: Comparison of the local moment calculated using 
FLEX and the ansatz-based technique presented in this pa- 
per against the exact Bethe-ansatz result. All schemes have 
the same weak-coupling behavior. Although FLEX overesti- 
mates the moment at strong coupling, the solution from our 
ansatz shows promising results. Both the underestimation of 
the moment that results from small cluster simulations, and 
the overestimation inherent in the FLEX approximation are 
corrected. 

comparison of certain quantities is possible. 

While the FLEX approximation predicts the correct 
AFM transition temperature (the Mermin-Wagner the- 
orem requires that there is only a transition at absolute 
zero), it is well known that FLEX describes local mo- 
ment formation incorrectly for couplings of the order of 
the bandwidth. The fundamental definition of the local 
moment is, 

3 

< fi >= 5(5+1) < (?i T -n x ) 2 >= -(< n > -2 < D >) 

where < D >=< n^n\ > is the expectation value of the 
double occupancy. 14] Since the potential energy term 
of the Hubbard model is UD, < D > can be extracted 
from the expectation value of the potential energy via, 

< V >= Tr[E(K, iu n )G(K, iw n )] = U <D> (9) 

Figure shows the local moment versus coupling, cal- 
culated using several different techniques (FLEX, the 
ansatz-based FLEX/QMC technique presented in this 
paper, and the exact ground state solution). There is a 
difference in temperature between the exact Bethe-ansatz 
solution of Lieb and Wu (T — 0) 2] and the approximate 
solutions computed here (T = 0.16VF) The temper- 
ature is still low enough to see the non-trivial effects of 
spatial fluctuations. The hopping, t = 0.25, fixes the 
the bandwidth of the non-interacting problem to unity 
(therefore the interaction energy scale for the boundary 
between weak and strong coupling regimes is U = 1). 
Calculations are carried out for values of U up to the 
bandwidth to examine behavior outside the perturbative 
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FIG. 6: Self energy in ID for various couplings. T = 0.08, 
N c = 4, N' c = 32, D = 1, n = 0.5. Each set of curves 
represents the variation of the self-energy across the Brillouin 
zone. Note how the momentum dependence becomes more 
pronounced as coupling increases. No causality problems are 
found for couplings of similar size to the bandwidth. 



regime. The exact result for a ID problem calculated 
from the Bethe ansatz solution is shown, and compared 
with the 4 site QMC, the FLEX solution for a cluster size 
of N c = 32, and the hybrid FLEX/QMC scheme. 

The results in figure |S] shows that all schemes have the 
same weak coupling behavior. The gradient of the low 
U curve increases with larger cluster size and decreased 
temperature to converge on the Bethe-ansatz solution. 
It can be seen that FLEX overestimates the moment at 
strong coupling, i.e. it underestimates the double oc- 
cupancy. Alternatively, results from the 4 site QMC 
approximation are inclined to underestimate the local 
moment at intermediate coupling, probably because the 
mean-field nature of the Nc = 4 DCA used here pre- 
dicts a metal insulator transition at a critical coupling 
U = Uc, rather than at U = + as expected in ID. The 
solution from the new ansatz shows promising results. 
The moment has the same weak coupling behavior seen 
in all the presented approximations. More importantly, 
the overestimation of the moment that was predicted by 
the FLEX approximation has been corrected. Remark- 
ably, the new hybrid method predicts a moment that 
closely follows the exact solution, and the mean- field like 
behavior of the DCA is greatly reduced. 

The self energy predicted by the hybrid scheme is 
shown in figure Calculations are also carried out for 
the ID system for a series of couplings. In this case, the 
temperature is T — 0.08W. Groups of curves with sim- 
ilar asymptotic behavior belong to the same coupling, 
and differ only by their location in the Brillouin zone. 
As coupling increases, the self-energy gets larger, as ex- 
pected. As compared to a conventional QMC calculation 
with Nc — 4, the self-energy has significantly more de- 
tail, representing many more points in the Brillouin zone. 





-0.01 




-0.02 - 




-0.03 




-0.04 


3 


-0.05 - 


£f 


-0.06 - 


E 






-0.07 




-0.08 ■ 




-0.09 ■ 




-0.1 - 




6 



This momentum dependence is most important for inter- 
mediate coupling results. 

By examining the U = 0.8W result, the advantage 
of the hybrid method over the traditional approach can 
be seen. With a QMC calculation, a linear interpola- 
tion would be used to extract the lattice self energy. 
This would mean that all the self-energy curves would 
be equidistant. In contract, the lines at the extremes 
(which correspond to the center and edge of the Bril- 
louin zone) are clearly much closer together. In fact, for 
QMC cluster sizes up to Nq = 6 in ID, only extremal 
behavior can be predicted, and therefore it is apparent 
that the Hybrid scheme generates results which are much 
closer to the true lattice self-energy. 



V. 2D MODEL 

An important quantity that can shed light on the 
physics of correlated electron systems is the existence 
and shape of the Fermi-surface. The Fermi-surface may 
be used as input for a variety of other theoretical tech- 
niques, and may be also be compared directly with re- 
sults from experimental methods such as angle-resolved 
photo-emission spectroscopy (ARPES) and de Haas van 
Alphen measurements. In this section, we examine fea- 
tures of the Fermi surface in 2D, by investigating the 
momentum-dependent electron density across the Bril- 
louin zone. The form of this quantity is closely connected 
to the shape of the Fermi-surface. 

The momentum dependent electron occupation may be 
calculated using the following formula, 



n k = G(k, iu> n ) exp(iw„0 + ) 



(10) 



where the lattice Green function has been constructed 
from the linear interpolation of the large cluster self- 
energy. The magnitude of the gradient of this quantity 
Vrik| is related to the Fermi-surface, since the gradient 
is largest at the Fermi-surface. 

The calculations in this section are performed at T = 
0.04 with coupling U = 2.0 = W. An N c = 4 QMC 
cluster is used throughout. 48 time slices were used for 
the Trotter decomposition of the QMC cluster. The 2D 
model has a simple cubic tight-binding dispersion with 
t x = t y = 0.25 and a small interplane hopping t z = 0.005 
to stabilize the solution. The Fermi-surface has been 
computed using the hybrid scheme for a series of fillings 
and large cluster sizes. 

One of the underlying aims of the hybrid scheme is 
to predict features consistent with very large systems, 
without the need for expensive QMC simulations. We 
first investigate the variation in |Vnk| as the cluster size 
is increased from N' c — 4 to N' c — 256, when the chem- 
ical potential is set to fi = 0.075 corresponding to a 2D 
Hubbard model just away from half-filling. The results 
of these calculations are shown in figure[7| As cluster size 
is increased, |Vrik| becomes sharper and better defined. 
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FIG. 7: Variation of the Fermi-surface, |Vn(k)| with cluster 
size. Calculations are performed with U = 2, T = 0.04, [i — 
0.075 and (a) iV^ = 4, (b) N' c = 16 and (c) N' c = 256. 
Larger cluster sizes lead to a more well defined surface, with 
a more electron-like character (see c). Contours are spaced 
every A|Vn(k)| = 0.05. 



Also a small number of finer features emerges. Although 
N' c = 64 is not shown, we note that the N' c = 64 cluster is 
close to convergence, with only small differences between 
N' c = 64 and N' c — 256. The larger cluster seems signif- 
icantly more electron- than hole-like, indicating that the 
band-gap inherent at half filling has become slightly nar- 
rower, although the diffuse nature of the surface indicates 
that the system is not in a Fermi-liquid state (this is to 
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FIG. 8: Evolution of the Fermi-surface in 2D, |Vn(k)|, for var- 
ious fillings calculated using the hybrid FLEX/QMC scheme. 
All calculations are performed with (7 = 2 = W, T = 0.04, 
N c = 4 and Nc = 256. Values of the chemical potential are 
(a) fj, = 0.025, (b) fi = 0.1 and (c) fx = 0.25. For values 
of filling close to half-filling, no Fermi-surface is seen, with 
a very diffuse gradient to the momentum dependent electron 
density, (b) shows an interim value, where the Fermi-surface 
begins to form as the edge of the gap is reached. It can be 
seen that the central peak develops greater intensity, and the 
beginnings of a curved Fermi-surface can be seen close to the 
(7r/2, 7t/2) point. Finally in (c) a clearly denned Fermi sur- 
face is seen, with significant additional intensity. Contours 
are spaced every A|Vn(k)| = 0.05 
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FIG. 9: Evolution of the filling, < n >, with chemical poten- 
tial, \x. The signature of a gapped state can be seen at low 
chemical potentials, with no variation of filling as chemical 
potential is increased. A gap of at least 0.05VF can be seen 
(W — 2 in this case). The three larger points correspond to 
the chemical potentials chosen in figure|H] The line is a guide 
to the eye. 

be expected close to half- filling, since there are significant 
spin- fluctuations) . 

An aspect of the 2D Hubbard model that is of general 
interest is the evolution of the Fermi-surface as filling is 
changed. In a dilute system, the Hubbard model is well 
described by the T-matrix approximation, which predicts 
Fermi-liquid behavior. At the other extreme, the half- 
filled system has a metal-insulator transition, with strong 
spin-fluctuations due to the proximity of a phase tran- 
sition to the anti-ferromagnetic state at absolute zero. 
Also, non-Fermi-liquid pseudogapped behavior has been 
reported for the 2D Hubbard model a little off half filling. 

In figure El |Vn k | is shown for the N c = 4, N' c = 256 
scheme. The results are computed for a fixed U and 
temperature value, and only the filling is varied. Rep- 
resentative graphs are shown to demonstrate the gapped 
state, the movement between gapped and metallic states, 
and the beginnings of the Fermi-surface. 

When close to half of the electronic states are occu- 
pied < n >— 1.0 (see the result for \i — 0.025), |Vrik| 
is quite disperse, with available electronic states across 
the Brillouin zone, and no clearly defined maximum as- 
sociated with a Fermi-surface. As filling is increased to 
< n >= 1.002 (/j = 0.1), the Fermi-surface becomes bet- 
ter defined, and some curvature due to free electrons is 
evident. For larger fillings, < n >= 1.02 (// = 0.25), the 
beginnings of an electronic Fermi-surface are evident. In 
addition to the electron-like states, an unusual hole-like 
surface can be seen. This is smaller in magnitude than 
the maximum corresponding to the electronic states, but 
it persists into the metallic state. 

In order to quantify the transition between gapped and 
metallic states, the expectation value of the occupation, 
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< n > , has been calculated as a function of the chemical 
potential, fi. The results are shown in figure |5J For 
chemical potentials of up to the order of /i = 0.1, the 
filling does not vary, sitting at < n >= 1.0 or half-filling 
(there is a slight error due to the QMC algorithm). This 
indicates the presence of a gap which is approximately 5% 
of the band width. Once the metallic regime is reached, 
the filling doesn't vary very quickly, which suggests that 
there are very few states close to the Fermi-energy. This 
may be an indication of pseudogapped behavior. 

VI. CONCLUSIONS 

We have introduced a new technique for the simula- 
tion of the Hubbard model. In this technique, short 
length-scale fluctuations are treated using QMC tech- 
niques, and the solution is supplemented with long length 
scale physics from the FLEX approximation. 

We have demonstrated that the new technique is suc- 
cessful at simulating the Hubbard model in ID and 2D, 
provided that the QMC cluster is sufficiently large. In 
particular, we find that results for the ID Hubbard model 
at low temperatures are very close to the exact Bethe 
ansatz solution. The new hybrid method works well for 
couplings up to the order of the band-width. This is im- 
portant, since spatial fluctuations on all length scales are 
expected to make a major contribution to the physics of 
the ID model. It is therefore expected that the hybrid 



QMC /FLEX calculations will give important insight into 
the ID problem, while remaining numerically cheap. In 
2D, the ansatz has been applied to the calculation of the 
Fermi-surface. We demonstrate that very large clusters 
of N' c — 256 can be simulated. For coupling of the order 
of the band- width, we demonstrate the evolution from 
insulating to metallic behavior. 

The new technique is important from two view points. 
First, QMC simulation of large clusters is very expensive, 
with computing time growing quickly with cluster size. 
Therefore, any improvement on the convergence proper- 
ties of the DCA method will aid in the accurate simula- 
tion of lattice models. Second, it has been proposed that 
a similar technique be used to combine DMFT with the 
ab-initio GW approximation |12|. In fact the two cluster 
solver approach is quite general, and it is expected that 
other models could be accurately treated using similar 
techniques. 
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